library(psych)
library(ggplot2)
library(GGally)
heart_data <- read.csv("/home/qgelena/code/scripts/R/stepik_projects/heart.csv")
dim(heart_data)
## [1] 303 14
names(heart_data)
## [1] "age" "sex" "cp" "trestbps" "chol" "fbs"
## [7] "restecg" "thalach" "exang" "oldpeak" "slope" "ca"
## [13] "thal" "num"
head(heart_data)
## age sex cp trestbps chol fbs restecg thalach exang oldpeak slope ca thal
## 1 63 1 1 145 233 1 2 150 0 2.3 3 0 6
## 2 67 1 4 160 286 0 2 108 1 1.5 2 3 3
## 3 67 1 4 120 229 0 2 129 1 2.6 2 2 7
## 4 37 1 3 130 250 0 0 187 0 3.5 3 0 3
## 5 41 0 2 130 204 0 2 172 0 1.4 1 0 3
## 6 56 1 2 120 236 0 0 178 0 0.8 1 0 3
## num
## 1 0
## 2 2
## 3 1
## 4 0
## 5 0
## 6 0
colSums(is.na(heart_data))
## age sex cp trestbps chol fbs restecg thalach
## 0 0 0 0 0 0 0 0
## exang oldpeak slope ca thal num
## 0 0 0 0 0 0
Я решила удалить строки с пропущенными значениями, так как данных для анализа достаточно
heart_data$ca[heart_data$ca == "?"] <- NA
heart_data$thal[heart_data$thal == "?"] <- NA
colSums(is.na(heart_data))
## age sex cp trestbps chol fbs restecg thalach
## 0 0 0 0 0 0 0 0
## exang oldpeak slope ca thal num
## 0 0 0 4 2 0
heart_data <- na.omit(heart_data)
Я решила изменить тип некоторых переменных и изменить их levels, так как в даном дата фрейме много векторов с качественными признаками
Качественные:
heart_data$sex <- factor(heart_data$sex)
levels(heart_data$sex) <- c("female", "male")
heart_data$cp <- factor(heart_data$cp)
levels(heart_data$cp) <- c("typical","atypical","non-anginal","asymptomatic")
heart_data$fbs <- factor(heart_data$fbs)
levels(heart_data$fbs) <- c("false", "true")
heart_data$restecg <- factor(heart_data$restecg)
levels(heart_data$restecg) <- c("normal","stt","hypertrophy")
heart_data$exang <- factor(heart_data$exang)
levels(heart_data$exang) <- c("no","yes")
heart_data$slope <- factor(heart_data$slope)
levels(heart_data$slope) <- c("upsloping","flat","downsloping")
heart_data$thal <- factor(heart_data$thal)
levels(heart_data$thal) <- c("normal","fixed","reversable")
heart_data$num <- ifelse(heart_data$num %in% c(0), 0, 1)
heart_data$num <- factor(heart_data$num)
levels(heart_data$num) <- c("Здоров", "Болен")
Дискретные: ca (не может количество сосудов быть дробным):
heart_data$ca <- factor(heart_data$ca) # не делаем конвертацию, так как это не обязательно;
Лайфхак: если сомневаетесь, дискретный признак или нет, постройте его гистограмму с расзмером бина = 1. Если признак непрерывный, то одинаковые значения повторяются очень редко. Если дискретный, то вы увидите несколько высоких столбиков.
Я построила гистограммы для всех количественных признаков с размером бина 1, судя по выше сказаному, эти признаки непрерывные.
ggplot(heart_data, aes(x = heart_data$thalach))+geom_histogram(binwidth=1)
ggplot(heart_data, aes(x = heart_data$age))+geom_histogram(binwidth = 1)
ggplot(heart_data, aes(x = heart_data$trestbps))+geom_histogram(binwidth = 1)
ggplot(heart_data, aes(x = heart_data$chol))+geom_histogram(binwidth = 1)
ggplot(heart_data, aes(x = heart_data$oldpeak))+geom_histogram(binwidth = 1)
Я использовала функцию describeBy, чтобы посмотреть на общие статистики, не учитывая данных с качественными признаками.
Можно заметить, что в выборке больше испытуемых, которые здоровы, 160 и 137, соответственно. В здоровых испытуемых можно увидеть крайние значения по переменной chol, сывороточный холестерин в отметке 564, что значительно больше, чем крайняя отметка для больных людей (409.0).
max кров. давление (thalach) в здоровых людей по среднему и медиане привышают значения в больных людей.
Также значения переменной oldpeak в здоровых людей ниже по медиане и среднему значению, чем в больных людей.
describe(heart_data)
## vars n mean sd median trimmed mad min max range skew
## age 1 297 54.54 9.05 56.0 54.73 8.90 29 77.0 48.0 -0.22
## sex* 2 297 1.68 0.47 2.0 1.72 0.00 1 2.0 1.0 -0.75
## cp* 3 297 3.16 0.96 3.0 3.29 1.48 1 4.0 3.0 -0.84
## trestbps 4 297 131.69 17.76 130.0 130.48 14.83 94 200.0 106.0 0.69
## chol 5 297 247.35 52.00 243.0 244.70 47.44 126 564.0 438.0 1.11
## fbs* 6 297 1.14 0.35 1.0 1.06 0.00 1 2.0 1.0 2.01
## restecg* 7 297 2.00 0.99 2.0 2.00 1.48 1 3.0 2.0 0.01
## thalach 8 297 149.60 22.94 153.0 150.91 22.24 71 202.0 131.0 -0.53
## exang* 9 297 1.33 0.47 1.0 1.28 0.00 1 2.0 1.0 0.74
## oldpeak 10 297 1.06 1.17 0.8 0.88 1.19 0 6.2 6.2 1.23
## slope* 11 297 1.60 0.62 2.0 1.54 1.48 1 3.0 2.0 0.51
## ca* 12 297 1.68 0.94 1.0 1.51 0.00 1 4.0 3.0 1.17
## thal* 13 297 1.84 0.96 1.0 1.79 0.00 1 3.0 2.0 0.33
## num* 14 297 1.46 0.50 1.0 1.45 0.00 1 2.0 1.0 0.15
## kurtosis se
## age -0.55 0.53
## sex* -1.44 0.03
## cp* -0.44 0.06
## trestbps 0.76 1.03
## chol 4.30 3.02
## fbs* 2.04 0.02
## restecg* -1.99 0.06
## thalach -0.09 1.33
## exang* -1.46 0.03
## oldpeak 1.44 0.07
## slope* -0.65 0.04
## ca* 0.19 0.05
## thal* -1.83 0.06
## num* -1.98 0.03
describeBy(heart_data[, -c(2, 3, 6, 7, 9, 11, 12, 13, 14)], heart_data$num, mat=T)
## item group1 vars n mean sd median trimmed
## age1 1 Здоров 1 160 52.643750 9.5511512 52.0 52.4531250
## age2 2 Болен 1 137 56.759124 7.8996701 58.0 57.2432432
## trestbps1 3 Здоров 2 160 129.175000 16.3739900 130.0 128.5390625
## trestbps2 4 Болен 2 137 134.635036 18.8967302 130.0 133.1171171
## chol1 5 Здоров 3 160 243.493750 53.7575499 235.5 239.2421875
## chol2 6 Болен 3 137 251.854015 49.6799374 253.0 251.0000000
## thalach1 7 Здоров 4 160 158.581250 19.0433045 161.0 160.0468750
## thalach2 8 Болен 4 137 139.109489 22.7106735 142.0 139.7837838
## oldpeak1 9 Здоров 5 160 0.598750 0.7871601 0.2 0.4617188
## oldpeak2 10 Болен 5 137 1.589051 1.3050061 1.4 1.4747748
## mad min max range skew kurtosis se
## age1 10.37820 29 76.0 47.0 0.1191477 -0.6652838 0.75508480
## age2 5.93040 35 77.0 42.0 -0.5736352 0.1352938 0.67491436
## trestbps1 14.82600 94 180.0 86.0 0.4346036 0.2568030 1.29447757
## trestbps2 14.82600 100 200.0 100.0 0.8014577 0.6514438 1.61445661
## chol1 44.47800 126 564.0 438.0 1.7032585 7.2146181 4.24990748
## chol2 50.40840 131 409.0 278.0 0.2771550 0.2207174 4.24444349
## thalach1 16.30860 96 202.0 106.0 -0.6789875 0.4109986 1.50550541
## thalach2 25.20420 71 195.0 124.0 -0.2857344 -0.3011150 1.94030378
## oldpeak1 0.29652 0 4.2 4.2 1.5772389 2.8336782 0.06223047
## oldpeak2 1.48260 0 6.2 6.2 0.7130820 0.2774105 0.11149420
Строим scatterplot для всех пар
pairs(heart_data[, -c(2, 3, 6, 7, 9, 11, 12, 13, 14)], col=heart_data$num)
ggpairs(heart_data[, -c(2, 3, 6, 7, 9, 11, 12, 13, 14)])
#Но в таком виде совершенно ничего не видно :(
#Советую узнать, как изменить размер графика и сделать его больше
ggpairs(heart_data, lower = list(combo = wrap("facethist", binwidth = 1, mapping=aes(colour=num, alpha=0.4))))
Было замечено несколько аутлаеров, некоторые из них при этом находятся в категории “здоровые”. Например, женщина с chol больше 500 и thalach около 160. Я бы уточнила у авторов исследования, не являются ли подобные значения ошибками измерений
pairs(heart_data[ ,c(1,4,5,8,10)], col = heart_data$num)
par (xpd=TRUE)
legend("bottomright", fill=unique(heart_data$num), legend=c(levels(heart_data$num)))
Хорошо заметна разница между группами здровых и больных по показателям oldpeak, thalach. У здоровых выше показатель thalach, у больных выше oldpeak
pairs(heart_data[ ,c(1,4,5,8,10)], col = heart_data$slope)
pairs(heart_data[ ,c(1,4,5,8,10)], col = heart_data$cp)
Можно заметить, что аутлаеры по переменной oldpeak лежат в одной (зелёной) категории slope. Также они лежат в одной (синей) категории cp
Разделения данных на чёткие облака отмечено не было
#“аутлайнеры по переменной oldpeak лежат в одной (зелёной) категории slope” а какому значению соответсвует зеленый slope?
pairs(heart_data[ ,c(1,4,5,8,10)], col = heart_data$slope)
par (xpd=TRUE)
legend("bottomright", fill=unique(heart_data$slope), legend=c(levels(heart_data$slope)))
гистограммы (не забудьте поэкспериментировать с размером бинов и объяснить, как изменился вид гистограммы) ящики с усами и с точечками вайлин плоты с точечками попробуйте facet_grid разбивайте графики по категориальным переменным и раскрашивайте их в разные цвета (попробуйте разные цветовые схемы)
ggpairs(heart_data[ ,c(1,4,5,8,10)], aes(col=heart_data$num, alpha = 0.4))
У здоровых выше показатель thalach, у больных выше oldpeak
ggpairs(heart_data[ ,c(1,4,5,8,10)], aes(col=heart_data$slope, alpha = 0.4))
“зелёный” slope - это “2” Value 2: flat" - плоский угол
ggpairs(heart_data, aes(col=heart_data$num, alpha = 0.4))
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
ggplot(heart_data, aes(x=num, y=oldpeak)) +
geom_boxplot(outlier.colour="red", outlier.shape=8,
outlier.size=4) +
geom_jitter(shape=16, position=position_jitter(0.2))
На этом графике мы видим, что по показателю oldpeak у больных, например, существенно больше медиана, а также “намного больше сам ящик”, то есть выше среднеквдратическое отклонение (Надеюсь, я не путаюсь в том, что ширина box plot позволеят судить о среднеквадратическом отклонении, раз для его построения мы откладываем 1,5 квартиля, размер которых как раз зависит от sd)
#{r} #ggplot(heart_data, aes(x=age, y=oldpeak, col=thalach)) + geom_point() + facet_grid(vars(num), vars(sex)) #
На этом графике видно, что больных женщин в испытании было существенно меньше. Если игнорировать outliers, то можно сказать, что характер распределения признака среди больных и здоровых похож у обоих полов
ggplot(heart_data, aes(x=num, y=oldpeak, col=sex)) +
geom_violin() +
geom_jitter(shape=16, position=position_jitter(0.2))
Различия в распределении признака идут по параметру здоровый/больной, а не полу, что различаются медианы и отклонение
ggplot(heart_data, aes(thalach, fill = sex, alpha = 0.5)) +
geom_histogram(binwidth = 5)
На этой гистограмме мы видим, что у женщин показатель thalach в среднем выше
ggplot(heart_data, aes(thalach, fill = sex, alpha = 0.5)) +
geom_histogram(binwidth = 50)
Давайте поэкспериментируем с bin`ом. Можно заметить, что при неправильно подобранном шаге информативность гистограммы существенно снижается
#я бы еще добавила легенду к цветам, чтобы было еще понятнее
ggplot(data=heart_data[!is.na(heart_data$thalach),], aes(num, thalach)) +
geom_violin(aes(fill = num), trim = FALSE, alpha = 0.3) +
geom_boxplot(aes(fill = num), width = 0.2, outlier.colour = NA) +
theme(legend.position = "NA") + geom_jitter(shape=16, position=position_jitter(0.2))
heart_data$new.var <- heart_data$age > 50
ggplot(heart_data, aes(x=chol, y=slope)) +
geom_point(alpha=0.5) +
facet_grid(new.var ~sex) +
scale_color_gradient(low = "blue", high = "red") +
theme_bw()
ggplot(heart_data, aes(x = chol)) + geom_histogram(fill='blue', col="yellow", binwidth = 7)
hist(heart_data$chol)
# Так никогда не делать (!):
# ggplot(heart_data, aes(x=oldpeak, fill =slope)) + geom_dotplot()
ggplot(heart_data, aes(x=num, y=chol)) + geom_boxplot()
В чанках, которые кидают ворнинги используй waning =F и message = F, чтобы не захламлять html. Например, вот тут стоило отключить ворнинги {r} library(psych) library(ggplot2) library(GGally)
Про количественные переменные нужно было понять дискретные они или непрерывные.
Вот эта операция совсем некорректная heart_data_dis <- heart_data[heart_data$num == "0" | heart_data$num == "1", ] 0 кодирует здоровых людей, а 1,2,3,4 больных. Тебе нужно было склеить всех больных людей, а ты вместо этого выкинула тех, у кого heart_data\(num == "2" или heart_data\)num == “3” или heart_data$num == “4”, то есть выкинула большой кусок датасета.
Там, где рисуешь парные скаттерплоты нужно описывать, что ты на них видишь – зависимости, неоднородности в данных, аутлаеры, …
Здорово, что научилась строить такой график ggpairs(heart_data_dis, lower = list(combo = wrap(“facethist”, binwidth = 1, mapping=aes(colour=num, alpha=0.4)))) Но в таком виде совершенно ничего не видно :( Советую узнать, как изменить размер графика и сделать его больше
Вот этот график классный! ggplot(data=heart_data_dis[!is.na(heart_data_dis$thalach),], aes(num, thalach)) + geom_violin(aes(fill = num), trim = FALSE, alpha = 0.3) + geom_boxplot(aes(fill = num), width = 0.2, outlier.colour = NA) + theme(legend.position = “NA”) + geom_jitter(shape=16, position=position_jitter(0.2)) я бы еще добавила легенду к цветам, чтобы было еще понятнее